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ABSTRACT 

The inclusion of rotationally-induced mixing in stellar evolution can alter the structure and compo- 
sition of presupernova stars. We survey the effects of progenitor rotation on nuclcosynthctic yields in 
Population III and II supernovae using the new adaptive mesh refinement (AMR) code CASTRO. We 
examine piston-driven spherical explosions in 15, 25 and 40 M Q stars at Z = and 1CP 4 Zq with 
three explosion energies and two rotation rates. Rotation in the Z = models resulted in primary 
nitrogen production and a stronger hydrogen burning shell which led all models to die as red super- 
giants (in contrast to the blue supcrgiant progenitors made without rotation). On the other hand, the 
Z = 10~ 4 Zq models that included rotation ended their lives as compact blue stars. Because of their 
extended structure, the hydrodynamics favors more mixing and less fallback in the metal free stars 
than the Z = 10~ 4 models. As expected, higher energy explosions produce more enrichment and less 
fallback than do lower energy explosions, and at constant explosion energy, less massive stars produce 
more enrichment and leave behind smaller remnants than do more massive stars. We compare our 
nuclcosynthctic yields to the chemical abundances in the three most iron-poor stars yet found and 
reproduce the abundance pattern of one, HE 0557-4840, with a zero metallicity 15 M , 2.4 x 10 51 
erg supernova. A Salpeter IMF averaged integration of our yields for Z = models with explosion 
energies of 2.4 x 10 51 ergs or less is in good agreement with the abundances observed in larger samples 
of extremely metal-poor stars, provided 15 M Q stars are included. Since the abundance patterns of 
extremely metal-poor stars likely arise from a representative sample of progenitors, our yields suggest 
that 15-40 M Q core-collapse supernovae with moderate explosion energies contributed the bulk of the 
metals to the early universe. 

Subject headings: Supernovae, nucleosynthesis, first stars 



1. INTRODUCTION 



Early chemical enrichment of the universe began 

with the deaths of the first stars at z ~ 20 

30 (iMackev et a .1 120031: iKitavama fc Yoshidal 120051: 



Greif et all l2007t iKarlsson et all 120081: IWhalen et all 
20081 ). The nature of the first stars, and hence the 
primordial initial mass function (IMF), is yet to be 
observationally constrained, but numerical models that 
proceed from well-posed cosmological initial conditions 
sug gest that they are very massive, from 25 - 500 
M^ feromm et all 119991 120021 iNakamura fe Umemural 
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20011: lAbel et all boOOL [2002t lO'Shea fe Normanl 120071 
Wise fc Ab el 2007) . The non-rotat ing stellar evolution 
models of iHeger fc Wooslevi (|2002l ) predict that 15 - 
40 M Population III (Pop III) stars die in conven- 
tional core-collapse supernovae, while 140 - 260 M Q stars 
explode as much more energetic pair-instability super- 
novae. For this reason, 15 - 40 M Pop III stars are 
considered to be "low- mass" primordial stars and 140 - 
260 M Q stars are termed "very massive" stars. Although 
direct observation of the first stars is not yet possible, 
clues to their nature may be extracted from their nucle- 
osynthctic imprint on second-generation stars. The gas 
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enriched by Population III stars is thought to be able to 
form low-mass stars that may still exist, retaining the 
chemical imprint of their progenitors. Fossils of this sec- 
ond generation of stars are now being sought in surveys 
of extremely metal-poor (EMP) an d hyper metal-poor 
(HMP) stars in the G alactic halo (jBeers fc Christliebl 
[20051: iFrebel et~aT1[2005l) . HMP stars, with [Fe/H] < -4, 
are thought to be enriched by only one or a few SNe. 
EMP stars, with -4 <[Fe/H]< -3, are more plentiful 
and show a smaller scatter in their abundance ratios, im- 
plying that they formed from gas that had been enriched 
by a representative sample of SNe. 

Many attempts have been made to reproduce the 
abundance patterns in EMP and HMP stars by mod- 
eling the evolution and explosion of mas sive stars, and 
comparing the yields to observations dNomoto et al 



2006 ; IUmeda fc Nomoto 2002, 2003, 2005: Uoggerst et al 



20091 ). and in doing so, gain a better understanding of the 
nature of the first stars. Some studies have focused on 
the potential role of very massive stars that die as pair- 
instability supernovae; others suggest that very energetic 
hypernova explosions may be responsible for these abun- 
dance patterns. Although constraints have been placed 
on the energy and origin of the explosion within the star 
and theoretical progress has been made dBurrows et all 
120061 : iMurohv fc Burrowsl l20?)l iMarek fc Jankal I2009D . 
the explosion mechanism that operates in this first gen- 
eration of stars is not well understood. Furthermore, 
core-collapse SNe leave compact remnants onto which 
the outer layers of the star may fall, preventing their es- 
cape and enrichment of the next generation of stars. The 
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degree to which the layers of the star are mixed by the 
Ray leigh- Taylor (RT) instability prior to shock breakout 
from the surface of the star, together with the amount of 
material that falls onto the central remnant, governs the 
final yield of the SN. 

The large one -dimensional explosi o n surveys per- 
formed to date (llwamoto et al.l 120051 : iTominaga et al.l 
l2007t iHeeer fc Wooslevl 12008ft are forced to approxi- 
mate these inherently multidimensional processes by ar- 
tificially mixing the layers of the star and selecting a 
"mass cut" that denotes the line between the inner ma- 
terial that falls back onto the central compact object 
and the outer material that escapes. Two-dimensional 
models have been constructed to obtain more accurate 
chemical yields for core-collapse SNe and to constrain 
one dimensional parameterizations of these explosions. 
Tominaga ( 2009) investigated one explosion mechanism, 
jet-induced explosions, with two-dimensional simulations 
of a 40 Mq primordial star and found good agreement 
using angle-delimited yields to the abundance pattern in 
one HMP star, but was not as sucessful in reproducing 
the abunace patterns in EMP stars. 

Previous studies of mixing in core-collapse super- 
novae were largely focused on replicating observations of 
SN1987A. The most s ucess ful of these were the st udies of 
iKifonidis etldl (|2003t) and iKifonidis et al.l (|2006t) . which 
modeled the explosion in more detail and from earlier 
times than previous studies. 

The first paper in this series (jJoggerst et al.l l2009f ) 
modeled spherical explosions of non-rotating 15 Mq and 
25 M Population III stars with the FLASH code. These 
stars ended their lives as compact blue supergiants, and 
experienced far less mixing and more fallback than the 
larger-r adius solar metallicity red supergiants also mod- 
eled in iJoggerst et al.l (|2009l ) . The non-rotating Popu- 
lation III model yields achieved only partial agreement 
with the abundances in the most metal-poor stars ever 
detected. In particular, they did not reproduce the 
amount of nitrogen observed in HMP stars because it 
was not present in the initial pre-supernova progenitor 
models in sufficient quantities. The zero-metallicity ex- 
plosions in this suite also manifested a high degree of 
fallback and a low degree of mixing, less than what had 
previously been assumed in the one-dimensional models. 
These two factors resulted in even greater [C+O/Fe] ra- 
tios than those in the HMP stars, and failed to match 
the abundance patterns observed in the larger sample of 
EMP stars. 

Replicating abundances, particularly of N, in EMP 
and HMP stars may require rotationally-induced mix- 
ing prior to the destruction of the star. Recent simula- 
tions of Popul ation III stellar evolut ion up to the point 
of explosion by Ekstrom ct al. ( 2008) have demonstrated 
that rotationally-induced mixing boosts burning in the 
H-shell after the star has left the main sequence, increas- 
ing 14 N production by as much as a factor of 10 6 and ex- 
panding the outer envelope of the star. The larger outer 
envelope in rotating models will impact post-explosion 
dynamics in the star, allowing mixing to occur on longer 
timcscalcs and perhaps producing more enrichment than 
in non-rotating models. 

In this, the second of a series of papers on elemental 
yields and mixing in Population III supernovae, we ex- 
amine rotating models in two-dimensional simulations to 



determine if they are responsible for the elemental pat- 
terns discovered in the EMP/HMP star surveys. These 
studies may also be used to constrain estimations of 
mixing and fallback for large ensemb les of one dimen- 
sional supernova models like those of iHeger &: Wooslevl 
(|2008l ). While these rotating progenitor models do not 
take rotationally-induced mass loss into account, and our 
explosions are spherical and not asymmetric (as may 
well be the case) these effects may have less impact 
on the stellar yield than the radius of the star at the 
time of its death. Our rotating models die as signifi- 
ca ntly larger, redder st ars than the non-rotating models 
of IJoggerst et al.1 (|2009l ). and we anticipate that this will 
lead to significantly enhanced mixing and reduced fall- 
back, and thus higher yields, as compared to non-rotating 
models. 

This survey is more complete in other respects. Our 
models capture mixing and fallback for a wider variety 
of progenitor masses, explosion energies, rotation rates, 
and metallicities than the previous survey. We include 
odd-numbered elements that were not included in our 
first survey, allowing for more detailed comparisons with 
EMP and HMP abundances. The greater number of el- 
ements also lays the groundwork for future light curve 
and spectroscopic studies of these supernovae. Our cur- 
rent study employs a new code, CASTRO, which can better 
handle the larger number of elements and will be more 
easily adapted to future studies than FLASH. 

In §[2] we describe the new multi-dimensional Eulerian 
hydrodynamics code CASTRO used to model the explo- 
sions. We describe the one-dimensional KEPLER models 
on which our simulations are based and outline our pa- 
rameter survey in § [HI and analyze the chemical yields 
of our simulation campaign in § |4] These yields are dis- 
cussed and compared with chemical abundance measure- 
ments in HMP and EMP stars in §0 In §[6] we conclude. 

2. NUMERICAL ALGORITHMS 

The explosion models in this study were implemented 
in two stages. First, one-dimensional supernova profiles 
for the progenitors were computed in the KEPLER code 
to capture all explosive nucleosynthetic burning. These 
profiles were then mapped onto a two-dimensional R- 
Z (axisymmetric) grid in the new CASTRO hydrodynam- 
ics code and evolved out to radii where Ray leigh- Taylor 
(RT) mixing ceased and the star was expanding essen- 
tially homologously. We computed thirty-six such mod- 
els, covering 3 masses, 3 explosion energies, 2 metallici- 
ties, and 2 rotation rates. 

2.1. CASTRO 

CASTRO l|Almgren et al.ll2010l) is a multi-dimensional 
Eulerian AMR hydrodynamics code. Time integration 
of the hydrodynamics equations is based on a higher- 
order, unsplit Godunov scheme. CASTRO can perform 
calculations in one dimensional radial, two dimensional 
cylindrical, or three dimensional Cartesian coordinates; 
all simulations in this paper were performed in two di- 
mensions with cylindrical coordinates. 

2.1.1. Equation of State 

CASTRO can follow an arbitrary number of isotopes 
or elements. The atomic weights and amounts of these 
elements are used to calculate the mean molecular weight 
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of the gas required by the equation of state. We followed 
the elements from hydrogen through zinc, so that our 
elemental yields could be compared with observations 
of abundances in metal-poor stars. 56 Ni was followed 
separately in order to calculate the energy deposited by 
radioactive decay of 56 Ni and 56 Co. 

The equation of state in our simulations assumed com- 
plete ionization and included contributions from both ra- 
diation and ideal gas pressure: 

p = /( P ,T)i a r 4 + M> (1 ) 

3 rripfj, 

E = f( P ,T)^ + 1.5^, (2) 
P m P ii 

where P is the pressure, a is the radiation constant, ks 
is Boltzmann's constant, T is the temperature, p is the 
density, m p is the proton mass, p is the mean molecular 
weight, and E is the energy. The function f{p,T) is 
a measure of the contribution of radiation pressure to 
the equation of state. It is 1 in regions where radiation 
pressure is important, i. e. where gas is optically thick, 
and in regions where radiation pressure is unimportant, 
i. e. where gas is optically thin, with a smooth transition 
in between. The function f(p,T) takes the form 



1 



f(P,T) 



if p > 10 9 gm cm 3 
or T < T, 



neg 



f(T) 



if p < 10 9 gm cm 3 
and T > T, 



where T neg is the temperature at which contributions 
to the pressure from radiation are negligable at 100 times 
less than that contributed by ideal gas pressure: 

Tneg = 3 ^ (3) 
lOOnipfia 

Above this temperature, contributions to the pressure 
from radiation will begin to dominate the equation of 
state, even though radiation pressure is unimportant in 
optically thin regions. Damping the radiation component 
of the EOS in optically thin regions prevented it from 
dominating in regions where it is not important. 

2.1.2. Radioactive Decay of 56 Ni 

Energy from the radioactive decay of 56 Ni to 56 Fe was 
deposited locally at each mesh point. The energy depo- 
sition rate from the decay of 56 Ni to 56 Co was: 

dE™ m = A5 BNi X56 Ni e- A56 N. t , 7 ( 56 Ni). (4) 
The decay rate of 56 Ni, As6 Ni , is 1.315 x 10~ 6 s -1 , and 
the amount of energy released per gram of deca ying 56 Ni 
is q(^ 6 Ni), which we set to 2.96 x 10 16 erg g _1 (|Wooslevl 
119881 ). Xse Ni is the mass fraction of 56 Ni in the cell. The 
mass fraction of 56 Co at a given time can be expressed 
in terms of the mass fraction of initial 56 Ni by 

As6 Nj 



Xs 



58 Co 



-Xs6 Ni (e 



- ^56 r 



so that the energy deposition rate from 56 Co is 

j 771 A56 N 



1, (5) 



e- A56 co*))A 56Co(? ( 
(6) 

x 10~ 7 s" 1 and 

an energy per gra m of decaying 5 6 Co, g( 56 Co), equal to 
6.4 x 10 16 erg g" 1 (jWooslevlll98"8h . 



6 Ni 



We used a decay rate A5e Co = 1.042 x 10 s 



2.1.3. Gravity 

CASTRO supports several different approaches to solv- 
ing for self-gravity. In the calculations presented here, 
we applied the monopole approximation for gravity. 
A radial average of the density was taken from the 
two-dimensional grid to create a one-dimensional den- 
sity profile. This profile was then used to compute a 
one-dimensional gravitational potential which was then 
mapped back onto the two-dimensional grid. Since per- 
turbations from spherical symmetry in the densities are 
very small, this approximation has a negligible effect on 
the final state of the calculation. A comparison between 
the monopole gr avity solver and the fu ll multigrid solver 
can be found in lAlmgren et al.l (|2010t ) . The radial ap- 
proximation to the gravitational potential has the advan- 
tage of being nearly as accurate as a multigrid solution 
for these explosions, but is much faster. 

Gravity from a point mass located at the origin was 
also included in the gravitational potential. The point 
mass represents the compact remnant left behind by 
the SN explosion. As inf ailing matter crosses the zero- 
gradient inner boundary near the origin, it is added to 
this point mass. 

2.1.4. CASTRO AMR 

The adaptive mesh refinement algorithm in CASTRO 
uses a nested hierarchy of logically-rectangular grids with 
simultaneous refinement of the grids in both space and 
time. The integration algorithm on the grid hierarchy is 
a recursive procedure in which coarse grids are advanced 
in time, fine grids are advanced multiple steps to reach 
the same time as the coarse grids, and the data at differ- 
ent levels are then synchronized. 

There is a regridding step in which increasingly finer 
grids are recursively embedded in coarse grids until the 
solution is sufficiently resolved. CASTRO uses as its de- 
fault r efinement criteria the "error estimator" of lLohnerl 
(| 1987T) , which is essentially the ratio of the second deriva- 
tive to the first derivative at the point at which the er- 
ror is evaluated. Details of the impleme ntation in vari- 
ous ge ometries are discussed more fully in lAlmgren et al.l 
(|201Q[) . The result is a dimensionless, bounded estimator, 
which allows arbitrary variables to be used with preset 
error indicators. 

We used density, pressure, velocity, and the abun- 
dances of 56 Ni, He, and O as our refinement variables. 
We slightly modified the refinement criteria for elemen- 
tal abundance so that they were applied only in regions 
where the abundance was greater than 10 -3 . Regions 
with abundances of an element lower than 10 -3 were not 
marked for refinement on the basis of elemental abun- 
dance alone. 

We also modified the refinement criterion for density. 
In order to control extraneous refinement as the star ex- 
panded on the grid, regions below a certain density were 
not marked for refinement below a certain refinement 
level. We chose a minimum "mass" at which refinement 
could occur, which was the same for all the simulations. 
56 CoTJic minimum density as a function of refinement level 
at which a cell can be tagged for refinement is given by 

Pmin,i = 2 2l /A x 10 20 gm cm- 3 , (7) 

where p m in,i is the minimum density at which refinement 
can occur at a given refinement level i and A a is the area 



4 



TABLE 1 

Blast Models: Total Angular 
Momentum and Fraction of Critical 
Angular Velocity 



model 


Ltot (gm cm 2 s 1 ) 1 




zl5-5 


6 1407E+51 


15.5 


z25-5 


6.1407E+51 


6.2 


z40-5 


6.1407E+51 


2.8 


ul5-5 


6.1407E+51 


11.3 


u25-5 


6.1407E+51 


4.6 


u40-5 


6.1407E+51 


2.1 


zl5-10 


1.2281E+52 


31.5 


z25-10 


1.2281E+52 


12.4 


Z40-10 


1.2281E+52 


5.5 


ul.5-10 


1.2281E+52 


23.3 


u25-10 


1.2281E+52 


9.2 


u40-10 


1.2281E+52 


4.1 



"Total angular momentum with which 
model was intialized 

'Percent critical angular velocity of the 
star this amount of angular momentum rep- 
resents 



of the cell at the Oth refinement level. The factor of 10 20 
was chosen because it kept refinement at the outer edge 
of the star to a reasonable level, where it was resolved 
with approximately the same number of cells as the outer 
edge of the helium shell, without losing important detail 
in the RT instabilities. 

2.1.5. Initialization of Multidimensional Data 

In mapping the radial data from KEPLER onto the two- 
dimensional axisymmetric grid in CASTRO, special care 
was taken to properly resolve the key elements of the 
simulations: the shock, the elemental shells, and the Fe 
core. In particular, both the Fe core and the O shell 
were resolved with a minimum of 16 cells. Doubling the 
resolution produced no essential change to the solution. 
Mixing ceased at same time, and the final distribution of 
elements as a function of mass was the same for simula- 
tions in which the Fe core was resolved with 16 cell and 
32 cells. The one-dimensional mapping results in explo- 
sions that are spherically symmetric in CASTRO: low-order 
departures from spherical symmetry are not included in 
the initial conditions, and do not develop in the course 
of the simulation.. Our models therefore only capture 
asymmetries of mode greater than I = 1 or 2 in the ex- 
plosions. 

3. PROGENITOR MODELS 

To create the explosion profiles that were mapped into 
CASTRO, the on e-dimensional Lagrangian stellar evolution 
code KEPLER (jWeaver et al.l 119781: IWoosfev et all l2002t) 
was first used to evolve the progenitors through all sta- 
ble stages of nuclear burning, up until their iron cores 
became unstable to collapse. At this point SN explosions 
were artificially initiated by means of a piston at constant 
Lagrangian mass coordinate that moved though the star 
with a specified radial history. The models presented in 
this paper used a piston located at the radius where the 
entropy was equal to A.Oks /baryon, which corresponded 
roughly to the base of the oxygen shell. The models were 
then evolved until all nuclear burning was completed (~ 



100 s after the start of the explosion) but prior to the exit 
of the forward shock from the helium shell (and hence 
before the formation of the reverse shock). Energy gen- 
eration was followed with a 19-isotope network up to the 
point of oxygen depletion in the core of the star and 
with a 128 isotope quasi-equilibrium network thereafter. 
Rotation was included in the models, which had the ef- 
fect of increasing mixing in semiconvective regions of the 
star, i.e. regions that are stable by the Schwarzchild but 
unstable by the Ledoux criteria. 

Our supernova models were taken from the survey 
of Heger (in prep). We adopted progenitor masses 
of 15, 25, and 40 M , which are thought to be on 
the lower end of the mass scale for Population II I 
stars (e. g. lO'Shea fe Norman 2007: Y oshida et aTEool . 
These masses were chosen in part because they were 
considered in p revious studies of zero- metallicity stars 
(iJoggerst et alj 120091: iTominagal 120091: iTominaga et~aT1 
120071: iHeeer fe Wooslevll2008D . We examined both zero- 
metallicity and Z = 1Q~ 4 Zq SN progenitors. These two 
metallicitics correspond to Population III stars and to 
stars on the cusp of the transition to Population II, where 
the gas has enough metals to cool efficiently and form low 
mass stars that can survive to the c urrent epoch (e. g. 
iMackev et al.ll2003tlFrebel et al.ll2007f ). Important differ- 
ences in stellar evolution arise in these two metallicitics; 
in particular, the 10 -4 Zq stars are less affected by ro- 
tation and are actually bluer than their rotating Z = 
counterparts. 

We chose models with explosion energies of 0.6, 1.2, 
and 2.4 Bethe, where 1 Bethe = 10 51 ergs (1.2 Bethe 
is near the obse rved median e xplosion energy for core- 
collapse SNe iHamuvi |2003 |). The star s were as- 
sumed to have no mass loss. iNomoto et all (|2006h and 
lUmeda fe Nomotol (|2005l ) proposed that the first gener- 
ation of SNe may have been "hypernovae" on the basis 
of observed abundances in the most metal-poor stars. 
We instead decided to study models closer to the more 
commonly observed explosion energies in the nearby (ad- 
mittedly solar-metallicity) universe. 

Two fiducial rotation rates "R" were also considered: 5 
and 10% of the critical angular velocity of a 20 M© solar- 
metallicity star, or total angular momentum of 6.1407 x 
10 51 or 1.2281 x 10 52 gm cm 2 s _1 , respectively. Models 
were initialized with the same amount of angular mo- 
mentum at a given fiducial rotation rate. These rates 
were 2 - 30% of the critical velocity for the stars in our 
survey, depending on their individual mass and metal- 
licity. This percentage is listed for each progenitor in 
Table [TJ Although the percentage of the critical veloc- 
ity at a given "R" decreases with mass, this amount of 
angular momentum is sufficient to drastically change the 
presupernova structure of the model from a non-rotating 
counterpart. The initial angular momenta of our stars 
ensure that if they were of solar metallicity, they would 
fall at typical to sub-typical rotation rates. Solar metal- 
licity massive stars typically rotate at about 30% of crit- 
ical angular velocity. Rotation speeds for Population 
III and II stars are not currently known, so in the ab- 
sence of more compelling data it is reasonable to adopt 
values typical of the local universe. Rapidly rotating 
stars would be more likely to explode with a high de- 
gree of asymmetry, perhaps in the form of a jet, than 
more slowly rotating stars. The spherically symmetric 
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explosions we have constructed for this study are more 
consistent with the slowly rotating progenitors we chose 
for our sample, though even slowl y rotating progenitor s 



tor our sample, though even slowl y rotating progenitor s 
may explode as ymmetri cally (see Murphy et al.l (|2009t) ; 
iMarek fe Jankal (|2009t l: IScheck et al.1 fcoOfl ) and refer- 



ences therein). 

Rotation can dramatically change the structures of 
the progenitors, especially the ze ro-metallicity ones 
(jEkstrom et al.l 120081: iHirschil l2007t) . Rotationally in- 
duced mixing between the helium core and the base of 
the hydrogen shell introduces C, N, and O, which leads 
to CNO burning at the base of the hydrogen shell. This 
dramatically increases the rate of energy production, 
puffing up the outer layer of the star. It also changes 
the interior structure, resulting in a uniformly mixed 
helium-hydrogen layer outside the CO core for the zero- 
metallicity models. Including even a small amount of 
angular momentum in the models effectively turns com- 
pact blue zero-metallicity stars into red giants. This has 
a profound effect on their post-explosion hydrodynam- 
ics. Rotation produces less of an effect in stars with low 
metallicity 

The nomenclature for our models is as follows. Metal- 
licity is denoted by either a "u" or a "z" at the begin- 
ning of the model name: "u" signifies 10 -4 Zq while "z" 
means zero metallicity. Next is the mass of the star: 15, 
25, or 40 Mq. The explosion energy is denoted by cither 
an "B", "D", or "G", corresponding to 0.6, 1.2, and 2.4 
Bethe, respectively. The rotation rate R is designated 
after the explosion energy by either a 5 or 10, the two 
fiducial rates defined above. For example, Z25D-5 is the 
25 Mq, zero-metallicity model, exploded with 1.2 Bethe 
of energy and rotating at R= 5. 

4. RESULTS 

The final composition of the ejecta of a core-collapse 
SN is determined by the degree to which the layers sur- 
rounding the compact remnant are mixed by the RT in- 
stability and by the amount and composition of material 
that accretes onto the remnant. Mixing and fallback are 
determined by the presupernova structure of the star and 
by the details of its explosion. 

4.1. Initial Stellar Structure and the Onset of 
Rayleigh- Taylor Mixing 

At the end of its life, a massive star is composed 
of shells whose isotopic weight increases toward the Fe 
core in the center of the star. These layers become 
unstable to RT mixing during the explosion. Heav- 
ier material from deeper layers is mixed outward and 
lighter material from shallower layers is mixed inward. 
A reverse shock forms when th e outgoing shock encoun- 
ters a region of increasin g pr 3 (jHerant fe Wooslevlll994l : 
IWooslev fe Weaver! 1 1995|) . When the shock encounters 
densities that fall less steeply than r -3 , i.e. a region of 
increasing pr 3 , it decelerates. The deceleration of the 
forward shock reverses the direction of the pressure gra- 
dient, which slows down the layers interior to the shock as 
well. Shocked material piles up into a high density shell. 
The reverse shock forms at the inner boundary of the 
high-density shell of decelerated matter and propagates 
down into the star toward its ce nter, slowing down th e 
deeper, inner layers of the star (|Kifonidis et al.l [20031 ). 
The deceleration of the shock creates a steep pressure 



gradient in the opposite direction to the gravitational 
and density gradients. The pressure gradient may over- 
whelm the gravitational gradient, and in doing so trigger 
the formation of RT instabilities in the material. The 
RT instability grows until the reverse shock has passed 
by, at which point the material becomes stable again and 
the instabilities cease to grow exponentially. 

The presupernova structure of the star gives some in- 
dication of where the reverse shock, and hence the RT 
instability, is likely to form. Figure [T] shows pr 3 for the 
R=10 progenitors in this study, scaled to show the loca- 
tion and shape of the helium and oxygen shells as well. 
We will show the final mixed structures of these stars in 
section § 14.31 The structure of the progenitor for R=5 
is qualitatively similar. Shown for compa r ison is a 25 
M Z =Z e star from iHerant fe Wooslevi (fl99l . The 
Z =Z Q and Z=10~ 4 Zq models are similar in that pr 3 
does not increase dramatically until the hydrogen-helium 
shell boundary. The reverse shock that drives RT mix- 
ing will form at this boundary. The Z = models are 
different. In these models, the shell boost at the base of 
the hydrogen layer has led to so much convection that 
a distinct helium shell no longer exists, so although the 
Z = Zq models are as large in radial extent as the 
Z = Zq models, the former lacks the helium shell present 
in the latter. Instead, helium and hydrogen are nearly 
uniformly mixed beyond the oxygen shell. The quantity 
pr 3 increases sharply at the boundary between the oxy- 
gen shell and the helium-hydrogen shell, indicating that 
the reverse shock will form, and RT mixing will com- 
mence, at this boundary deeper within the star. 

4.2. Effects of Rotation on Mixing 

Figure [2] shows the abundance distribution of individ- 
ual elements as a function of mass for models zl5G-5 and 
zl5G-10. Although they are not identical, the distribu- 
tions are quite similar for the two models. The other 
R=10 SN models in our sample displayed a similar if 
not higher degree of congruence with their R=5 coun- 
terparts. If the star rotates, the degree of rotation has 
little effect on the final distribution of elements in the 
explosions, at least for the rotation rates explored in our 
survey. The reason for this is that once a Z — star 
acquires even a small amount of rotation and becomes a 
red g iant rather than a compact blue star (|Joggerst et al.l 
I2009T ) , the structure of the red giant is relatively insensi- 
tive to the degree of rotation. Rotating Z = stars are 
about an order of magnitude larger in radial extent than 
non-rotating Z = stars, though our 25 M rotating 
models are closer in size to their non-rotating counter- 
parts than 15 and 40 M Q models. Consequently, as the 
ejecta propagates through the star it encounters nearly 
the same gas profile whether it is an R=5 or an R=10 
star and becomes similarly mixed. Mixing within the R 
series in Z = 10~ 4 Zq stars is similar because the intro- 
duction of a modest amount of rotation has little effect on 
the structure of progenitor in these low metallicity mod- 
els. Because models of cither metallicity were so similar 
for R=5 and R=10, we will focus on the R=10 mod- 
els for the remainder of this paper. However, because 
the structures of red giants are so distinct from compact 
blue stars, mixing and fallback vary dramatically with 
metallicity, as we discuss next. 
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Fig. 1. — Pre-supernova structure for Z = (top row) and Z = 10 -4 (middle row) used in this survey, along with the presupernova 
structure of a solar-metallicity star (bottom row) for comparison. The helium and oxygen shells are shown in light grey and dark grey, 
respectively, and the quantity pr 3 is shown as a solid line. When the forward shock encounters a region of increasing pr 3 , a reverse shock 
forms, which leads to the formation of the R-T instability. For the z-series models, this region is the boundary between the helium-hydrogen 
and oxygen shells. For the Z = 10~ 4 Zq stars, the boundaries between both the helium and hydrogen shells and the oxygen and helium 
shells are unstable. 

gen shell, making the outer envelope of the star entirely 
convectivc as it left the main sequence. The outer radii 
of the Z = stars ranged from 6.7 x 10 12 to 1.32 x 10 14 
cm for the 15 to 40 M Q stars, respectively. The outer 
radii of the Z = 10~ 4 Z & stars were about an order of 
magnitude smaller, ranging from 2.4 x 10 12 cm for the 
25 M Q star to 4.18 x 10 12 for the 40 M star. Since the 
time it takes the reverse shock to propagate through the 
layers of the star sets the timescale on which RT mix- 
ing can occur, we expect that RT instabilities will have 
more time to develop, and greater mixing will occur, in 
the more "puffed up" Z — stars. Stronger explosions 
will induce more mixing than weaker ones. The strength 
of the reverse shock, the degree to which the pressure 
gradient is reversed, and therefore the violence of the RT 
instabilities depends on the strength of the initial shock 
and hence on the energy of the initial explosion. 

We evolved the initial KEPLER supernova models in 
CASTRO past the formation of the reverse shock and RT 
instabilities, until RT mixing ceased and the SN remnant 
was expanding nearly homologously. We found that the 
time at which RT mixing halted, and the degree to which 
it mixed together the isotopic shells of the model, indeed 
varies with mass, metallicity, and explosion energy. Re- 
sults for the 18 R=10 models are presented by mass in 
Figures [3JH and El 

4.3.1. 15 M Q Models 

In Figure [3] we show the mass abundance of helium, 
oxygen, silicon, and iron in 15 M stars of 0.6 (B, at 




2 6 10 14 2 6 10 14 
Mass (M Q ) 

Fig. 2. — Comparison of elemental abundance as a function of 
mass coordinate for star zl5G. The solid line is the model zl5G-5 
(see Table [TJ; the dotted line corresponds to model zl5G-10. No 
substantial difference in element distribution was found between 
the different rotation speeds for any stellar model. This model is 
representative of the average differences between stars at R=5 and 
R=10. 

4.3. RT Mixing in Models by Mass, Metallicity and 
Energy 

The radius of the outer envelope of stars of a given mass 
is quite different for the two metallicities. The Z = 
models are much redder than the Z = 10 -4 Zq mod- 
els because their outer envelopes are larger. This is due 
to the hydrogen shell boost that dramatically increased 
the rate of hydrogen burning at the base of the hydro- 
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Fig. 3. — Mass abundance of He, O, Si, and Fe in Z = (top) and 10 — 4 (bottom) 15 Mq stars after RT-driven mixing in the model 
has ceased. Z = stars, which die as large red giants, show much more mixing than Z = 10 -4 Zq stars, which die as smaller blue giants. 
The snapshots are of the simulation at 3.3 X 10 s s, 4.4 X 10 5 s, and 4.0 X 10 5 s for zl5B, zl5D, and zl5G, respectively, and similarly, at 
7.9 X 10 4 s, 6.1 X 10 4 s, and 5.0 X 10 4 s for ul5B, ul5D, and ul5G. Mixing increases with explosion energy, which is 0.6, 1.2, 2.4 Bethe from 
left to right across the panels. The jetting along the y- and x-axes, a numerical artifact, causes the departures from spherical symmetry in 
the u-series models but does not substantially affect the conclusions in this paper. 



left) 1.2 (D, center) and 2.4 (G, at right) Bethe explo- 
sions after the end of RT mixing. The snapshots are of 
the simulation at 3.3 x 10 5 s, 4.4 x 10 5 s, and 4.0 x 10 5 s 
for zl5B, zl5D, and zl5G, respectively, and similarly, at 
7.9 x 10 4 s,6.1 x 10 4 s, and 5.0 x 10 4 s for ul5B, ul5D, 
and ul5G. Beyond these times, the distribution of ele- 
ments as a function of mass does not change. The zero- 
metallicity models (top panels) show a higher degree of 
mixing than their Z = 10~ 4 counterparts, as expected. 
As explosion energy increases (from left to right in Fig- 
ure [3j) the degree to which the inner layers of the star 
are mixed also increases. This is most apparent in the 
u-series models. In model ul5B-10, the RT instability 
does not fully develop, leaving only a hint of "fingers" 
at the H-He and He-0 shell boundaries. In model ul5D- 
10, where the explosion energy has doubled, RT mixing 
has penetrated as far as the silicon shell, and silicon fin- 
gers extend well into the oxygen shell. Iron, however, 
remains unmixed, even in model ul5G-10. The Z = 
models show nearly complete mixing, however. Fe from 
the inner core is mixed out to at least as far as the oxygen 
shell, and the RT instability has become fully non-linear. 
The original fingerlikc morphology has been significantly 
blurred over the many mixing timescales over which the 
instability has had time to develop. 



4.3.2. 25 M Models 

We show mixing in the 25 Mq models in Figure |H 
Here, as at 15 M , RT mixing has had a longer time 
to develop in the z-series than in the u-series models. 
The instabilities have evolved to the point that the ini- 
tial fingers have been nearly erased in the Z series. The 
snapshots show the distribution of elements at 3.1 x 10 4 
s, 6.3 x 10 4 s, and 2.7 x 10 4 s for z25B, z25D, and z25G, 
and 1.4 x 10 4 s,5.3 x 10 4 s, and 1.2 x 10 5 s for models 
u25B, u25D, and u25G, respectively. The center of the 
star has become mixed, and iron has reached the outer 
layers of the star. The initial radii of the 25 Mq Popula- 
tion III and Population II models differ only by a factor 
of 3 or so, less than the difference in initial radii between 
15 or 40 Mq stars of different metallicity. This is why 
the amount of mixing seen in the 25 Mq models is more 
similar than it is for stars at 15 and 40 M Q . 

While the figures do not show it clearly because they 
have been cropped to show just the mixing in the star, 
"jetting" , or excessive flow along an axis, is present along 
the y-axis for the z-series models and is more pronounced 
in the higher explosion energy models (D and G) than 
in B. These jets extend as far as 2 x 10 14 cm, or twice 
the radius of the mixed region, for the z-series models. 
The u-series models again exhibit more mixing at the 
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Fig. 4. — Mass abundance of He, O, Si, and Fe in Z = (top) and 10~ 4 Zq (bottom) 25 Mg stars after the end of RT-driven mixing. 
The snapshots are of the simulation at 3.1 X 10 4 s, 6.3 X 10 4 s, and 2.7 X 10 4 s for z25B, z25D, and z25G, and 1.4 X 10 4 s, 5.3 X 10 4 s, and 
1.2 X 10 s for models u25B, u25D, and u25G, respectively. Red Z = stars again show much more mixing than blue Z = 10 — 4 Zq stars, 
although it is not as extreme as in the 15 Mq models, in which the difference in outer radius between the z- and u-series progenitors was 
greater. Mixing again rises with explosion energy, which is 0.6, 1.2, 2.4 Bethe from left to right across the panels. Spurious jetting is also 
visible along the y- and x-axes in the u-series models. Like the 15 Mg stars shown in Figure [3] both mixing and the amplitudes of the RT 
instabilities clearly increase with explosion energy at both metallicities. 



higher explosion energies, although the effect is not as 
pronounced as in the 15 M Q models. The deviation from 
spherical symmetry in the inner, unmixed shells of silicon 
and oxygen is a numerical artifact of the axisymmetric 
geometry, which is discussed in greater detail below in 
§ 14.3.61 The artifact becomes more prominent with ex- 
plosion energy. 

4.3.3. 40 M Q Models 

The 40 M Q mass models are shown in Figure [5] at 
4.3 x 10 5 s, 3.5 x 10 5 s, and 2.9 x 10 5 s for models 
z40B, z40D, and z40G, and 2.5 x 10 4 s, 2.6 x 10 4 s, and 
4.9 x 10 4 s for models u40B, 740D, and u40G, respec- 
tively. These explosions experience more fallback than 
the lower mass models, especially at explosion energies 
of 0.6 Bethe, where they leave behind remnants of 8.1 
and 9.1 M Q for z40B and u40B, respectively Mixing 
also occurs in less time. Even in the red giant Z = 
models, the Fe core has completely fallen back onto the 
remnant at the center of the simulation by the time RT 
mixing begins in all but the most energetic explosions, 
so these models produce virtually no iron. The variation 
in mixing with explosion energy is also more apparent 
in the 40 M u-series models than in the 15 and 25 M Q 
u-scries runs. 



4.3.4. Variation of Mixing with E ex , M, Z 

The final appearance of the SN ejecta at a given mass 
does not vary qualitatively with rotation rate, but mass, 
metallicity, and explosion energy do govern the final com- 
position of the ejecta and the mass of the remnant. Fig- 
ures (3J 21 and O give a good qualitative picture of the 
morphology of mixing in the explosion models; in Figure 
[6]we quantify the propagation of He. C, N, O, Mg. and Fe 
in the ejecta as a function of mass coordinate for the 18 
R=10 models after RT mixing has ceased. Several trends 
are apparent. Lower- mass models show less fallback and 
more mixing in the internal layers than higher-mass mod- 
els. The z-series SNe have far more mixing than u-series 
SNe. SNe with higher explosion energies exhibit more 
mixing and less fallback than SNe with lower explosion 
energies. In particular, the B series SNe with sub-normal 
explosion energies, 0.6 Bethe instead of the canonical 1.2 
Bethe, eject almost no iron with the exception of model 
zl5B. 

The z-series models all show more mixing than their 
u-series counterparts. The 25 M Q models show the most 
mixing of the models in the u-series, while the 40 M Q u- 
series runs show the smallest degree of mixing. All the 40 
M Q models experience a great deal of fallback, but the 
u-series models show the most because they are more 
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Fig. 5. — Mass abundance of He, O, Si, and Fe for Z = (top) and 1O _4 Z (bottom) 40 M stars at 4.3 X 10 s s, 3.5 X 10 5 s, and 
2.9 X 10 5 s for models z40B, z40D, and z40G, and 2.5 X 10 4 s. 2.6 X 10 4 s, and 4.9 X 10 4 s for models u40B, 740D, and u40G, respectively. 
The same trends in mixing with mass, metallicity, and explosion energy as in the 15 and 25 Mq stars are evident in these models. 

compact. The higher explosion energy models exhibit 
less fallback. 

4.3.5. Comparison with Kepler Estimations of Mixing 

The large one dimensional surveys of SNe derive fi- 
nal estimates of elemental yields by artificially mixing 
the layers of the SN after explosive nucleosynthesis is 
complete. Surveys employing the Kepler code estimate 
mixing by passing a running boxcar average of width 
(in mass coordinate) W through the star, where W is 
10% the mass of the helium core. That is, the abun- 
dances at points that fell within a bin of width W were 
averaged together and set to this average, the bin was 
moved forward by one point, and the process repeated, 
moving outward through the star. This is done 4 times, 
artificially mixing the mass shells. In Fig [7] we compare 
KEPLER estimations of mixing with our two dimensioanl 
CASTRO results. In our two-dimensional CASTRO simula- 
tions, we find that some elemental shells are more mixed 
than others. The RT instability typically forms at the 
He-H or O-He boundary and advances inward. This re- 
sults in the helium and oxygen layers being more mixed 
than in KEPLER and the iron, and sometimes silicon, lay- 
ers being less mixed than the KEPLER estimations for the 
Z-series models. Our compact U-series models show less 
mixing in all elements than in KEPLER. 

4.3.6. Numerical Artifacts & Model Limitations 



Numerical artifacts arising from the mesh geometry are 
most prominent in the higher explosion energy, u-serics 
models, but they are present in all the runs. Jetting is 
visible along both axes in the z-series models in Figures[3l 
HI and [51 but is more pronounced along the y-axis. This 
is a well-known artifact that appears on axisymmetric 
coordinate grids: PROME THEUS and FLASH (jFrvxell et al.l 
l!99lHCalder et al.ll200l ) also manifest it. The deviation 
from spherical symmetry (most clearly visible in models 
ul5D-10 and ul5G-10) is also a result of the tendency 
of the shock to be most exaggerated in the direction in 
which it is moving, this time with respect to the reverse 
shock. Ripples in the shells appear along both the x- and 
y-axes. 

These anomalies do not change the essential results 
of our simulations. The z-scrics models would show a 
greater degree of mixing even without the presence of a 
jet directed either outward or inward. In addition, be- 
cause of the axisymmetric geometry, the mass lying along 
the y-axis is only 1 /X (where X is the distance along the 
x-axis) times the mass at a similar radius from the ori- 
gin, so the jet contains little mass. Additional mixing in 
the form of material pushed out farther than it should be 
along this axis was removed from the averages simply by 
excluding a 5° cone of material around the axis. Fallback 
is not affected by jetting in our simulations. 

In addition to these minor numerical effects, there are 
some inherent limitations to our models. As discussed in 
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Fig. 6. — Distribution of elemental abundances as a function of mass. Stellar mass increases from left to right; the top three panels are 
Z = stars; the bottom three rows are Z = 10~ 4 Zq stars. Explosion energy increases from top to bottom for each mctallicity. Individual 
elements are coded by color. Zero-metallicity stars in general show more mixing than the slightly enriched stars. Higher explosion energy 
SNe exhibit more mixing and less fallback, and lower mass SNe have more mixing and less fallback. 
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Fig. 7.— Distribution of He, O, Si, and Fc for Z = (top) and 10~ 4 Z Q (bottom) 15 M (left), 25 M Q (center), and 40 M (right) 
for 1.2 X 10 51 erg explosions. CASTRO 2-d simulation results (solid lines) are compared with Kepler approximations of mixing (dashed 
lines). Kepler underestimates mixing of oxygen in the Z-series models, and overestimates mixing of elements heavier than oxygen. Kepler 
overstimates mixing for all elements as compared to the results of 2-d simulations of U-series models. 



11 



§ 14.3.61 these explosions are spherical, and capture initial 
perturbations ol only much higher order than I = 1 or 2. 
Jet-driven explosions, as well as those with asymmetries 
of order I = 1 or 2, may be common in the early universe, 
especially in rapidly rotating stars. Such expl osions may 
have very different nuclcosynthetic patterns (jTominagal 
l2009t) . By examining only more slowly rotating stars, it is 
likely that the mapping of one-dimensional profiles onto 
two-dimensional grids in CASTRO is not excluding serious 
low-order perturbations, but not enough is known about 
the explosion mechanisms of core-collapse SNe, let alone 
Population III and II SNe, to be certain. 

Simulations of spherical explosions in core collapse 
SNe, especially those of blue progenitors like the u- 
series Z = 10~ 4 Zq SNe presented here, have his- 
torically underestimated mixing and final elemen t al ve- 
locities dHerant fe Wooslevl 11994; iFrvxell et all 119911: 
iKifonidis et al.ll2003l l2006ft . Previous studies that at- 
tempted to match observations of SN 1987A (which died 
as a compact blue giant) with techniques similar to those 
in our study were largely unsucessful, finding less mix- 
ing and lower v elocities for the heavi er elements than 
were observed. IKifonidis et al.1 (j2006| ) came closest to 
reproduxing observations of SN1987A by adding addi- 
tional energy to the base of the explosion in two dimen- 
sions in the first few seconds after the explosion. We do 
not employ this technique in this paper. Solar metallic- 
ity models of red supergiants, like our Z = models, 
come closer to matching observations of similar red SNe 
(jHerant k Wooslevlll994j ). It is possible that more mix- 
ing would occur if the details of the first few seconds 
of the SN explosion were known with greater accuracy. 
Our mixing estimates should therefore be taken as lower 
bounds, though our mixing estimates for the red Z = 
Zq models in this paper are likely more accurate than 
for their blue Z = 10~ 4 Zq counterparts. 

5. DISCUSSION 

5.1. Comparison with Metal-Poor Stars 

Comparison of our theoretical yields to direct obser- 
vations of Z — and Z = 10~ 4 Zq SNe is not cur- 
rently possible, as these objects do not exist in the 
nearby universe. However, the fossilized yields of such 
explosions may be detectable in the abundance pat- 
terns of extremely metal deficient h alo stars. The thre e 
most iron-poo r stars, HE0557-4840 (jNorris et al.ll2007t ). 
HE0107-5240 (IChristlieb et al.l l200l. and HE1327-2326 
(jAoki et al.ll2006l ). are thought to have been enriched by 
at most one or a few SNe, while EMP halo stars with 
—4 <[Fc/H]< —3 are believed to be imprinted by ele- 
ments from a wider sample of early SNe. In Figure |8] we 
show our SN yields together with observed abundances in 
the three HMP stars and the average o f the abundances 
of t he larger samp l e of E MP stars in the lLai et al.l (|2008[ ) 
and ICavrel et al.l (|2004 ) surveys, which contain 15 and 
22 stars, respectively. 

We calculate our theoretical yields by assuming the 
ejected material is completely mixed, as is usual for this 
kind of calculation. All short-lived isotopes ( 56 Ni for 
example) are reported as their eventual decay products. 

Figure [8] shows that our SN models do not reproduce 
the high CNO abundances present in HE 01327-2326 and 
HE 0107-5240 when the ratio of Mg to Fe is held fixed. 
With its lower ratios of C and N to Fe, HE 0557-4840 



is a better match to some of our che mical yields, and 
is al so closer to the abu ndances in the lLai et all ()2008l ) 
and ICavrel et al.l (|2004f ) data. It is also clear from the 
diagram that some models cannot produce the [Fe/Mg] 
values of any of the observations. All but the most ener- 
getic of the 40 M Q SNe yield an [Fe/Mg] ratio too small 
to appear in Figure [5J The 25 M© explosion with ener- 
gies less than 2.4 Bethe show little overlap with either the 
HMP or EMP data, either. The u-series 15 M Q models 
produce a [Mg/Fe] ratio that is too high. 

5.1.1. HMP Stars 

While HE 0107-5240 and HE 1327-2326 show higher 
ratios of C,N, O and Na to Fe than our models, the yield 
from one model, zl5G, does fit the abundance pattern 
of halo star HE 0557-4840 with the exception of car- 
bon. Oxygen is an upper limit, and [C+O/Fe] in this 
star is far lower than in the other two HMP stars. The 
other two HMP stars, HE 0107-5240 and HE 1327-2326 
show much higher [C+0/Fe] and [C+0/Mg] ratios that 
our model yields do not replicate. We compare yields 
from the explosions in our survey that most closely re- 
produce [Mg/Fe] for these two HMP stars with their 
observed abundance ratios in Figure GO Model zl5G 
provides a very close fit to the abundance patterns in 
HE0557-4840. The only elements not in agreement are 
carbon, where our model underprcdicts [C /Mg] by about 
0.5 dex, and chromium, where our model underp r edicts 
[Cr/Mg] by a little more than 1 dex. iNorris et all (|2007f ) 
compare their observed abundances for this star to the- 
or etical yields from the "m ixing and fallback" model s 
of lUmeda fe Nomotol (|2003f ) and llwamoto et al.l (|2005l) . 
which overpredict [Cr/F e]. The hypernoya mod els of 
lUmeda fe Nomotol (|2005h and lNomoto et al.l (|2006| ) yield 
a [Co/Fe] ratio of 0.5, which is considerably higher than 
the observed value -0.3 for this star, or the value of -0.2 
predicted by model zl5G. 

We also note that the n o n-rota ting spherical explosion 
models in iJoggerst et al.l (|2009T> also fail to reproduce 
the abundance patterns found in metal-poor stars. In 
particular, far less nitrogen is formed in their models 
than observed or produced in the rotating models. This, 
not surprisingly, suggests that Population III stars were 
indeed rotating, but at the fairly low speeds modeled 
here the exact rate seems to be unimportant. 

Reconciling our models to the abundance patterns in 
HE0107-5240 and HE1327-2326 is more problematic. En- 
richment by a companion AGB star has not been ruled 
out for either star, and would result in abundances of C, 
N, O, and Na well above that with which the stars were 
initially born. This has been proposed a s an explanation 
for t he abundances in H E0107-5240 (IChristlieb et al.1 
2001 and HE132 7-2326 (jAoki et al.l |2006|), although 
Frebel et al.l (|2008D presents evidence that HE1327-2326 
is not a member of a binary system and thus could not 
have been enriched by an AG B companion. While the 
mixing and fallback models of llwamoto et al.1 (|2005f ) ac- 
count for the abundance patterns in these stars, they are 
still heavily tuned one-dimensional parameterizations, 
not multidimensional calculations with realistic physics. 
T he jet-driven expl osions of a 40 M Population III star 
bv lTominagal ([2009') do reproduce the features of the ear- 
lier mixing and fallback models, providing a mechanism 
for the mixing-fallback parameterization. Their angle av- 
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Fig. 8. — Comparison of elemental yield r atios from SN m odel s with those observ ed in the 3 most metal-poor stars and with average 
abundances for EMP stars from the surveys of Lai ct al. (2008) and Cayrel ct al. (2004). Observations are denoted by a star symbol; models 
are denoted by open symbols where the size of the symbol varies with explosion energy, shape with mass, and metallicity with color. None 
of our stellar yields fit the high [C+O/Fe] of HE1327-2326, nor are we able to closely reproduce the abundance pattern of HE0107-5240 for 
elements heavier than Na. The z!5D and z!5 yields arc good fits to HE0557-4840. 



eraged yields match the abundances observed in HE1327- 
2326 and HE0107-5240-in particular, they produce the 
high [CNO/Fc] ratios that our models do not. 

Our models never show [C/Mg] and [O/Mg] ratios 
greater than 0.25 and 0.5, respectively. HE0107-5240 
has [C/Mg] and [O/Mg] ratios of 2.4 and 1.4, respec- 
tively; the values of [C/Mg] and [O/Mg] for HE1327- 
2326 are 1.1 and 1.5. That our models underreproduce 
the combined [CO/Mg] values for these stars by at least 
2 dex does not imply that they are not representative of 
Population III SNe. Since only [C/H] and [O/H] ratios 
greater than -3.5 produce enough fine structure cooling to 
trigger th e rollover from Popu l ation III to II star f orma- 
tion (e. g. lMackev et~aT1l2003f ). iFrebel et al.1 (|2007| ) have 
argued that any star surviving until today with [Fe/H] 
< —4 must have enhanced abundances of C and/or O. 
Consequently, although perhaps much more plentiful, 
Population III SNe that did not produce high [C/Fe] 
and [O /Fe] ratios would not have imprinted their chemi- 
cal signatures upon a subsequent generation of stars that 
was long lived. More exotic explosion mechanisms may 
not have been predominant at high redshift, and may not 
even be necessary for reproducing the abundance pat- 
terns in the more iron-poor stars yet found. Fig [6] shows 
that for some stars, C, N, and O are mixed out to higher 
mass coordinate than Mg and heavier elements. It is 
possible that HE0107-5240 and HE1327-2326 could both 
have formed from gas that was preferentially enriched by 



these further-flung lighter elements. 

5.1.2. EMP Stars 

Many more stars have been found with [Fe/H] values 
between -4 and -3 than below -4.5. The abundances in 
these stars show less scatter than in the small sample of 
three below -4.5. Figure [5] shows the average ratios of 
various abund ances to Mg as a fun ction of [X/Mg] from 
the surveys of lLai et al.l (|2008l ) and lCavrel et all (|2004l ). 
Figure ITOl shows IMF averaged yields for 0.6, 1.2, and 2.4 
Bethe explosions and Z = and 10 -4 Zq. We used a 
Salpeter IMF, 

N(M)dM oc M~ a (8) 

and the cannonical value for a of 2.35. 

The fact that the u-series models don't produce enough 
nitrogen should perhaps be ignored, since these models 
were constructed with yields from non-rotating Z = 
stars, which did not produce nearly the nitrogen that was 
formed in the rotating Z = Z G models. Explosions of 
1.2 and 2.4 Bethe (i.e. "D" and "G" models) for Z = 
come closest to matching the abundances of the EMP 
stars. 

An IMF-average over yields from the Z = models 
provides a good fit to the Cayrel and Lai data. Alu- 
minum is over produced by about 0.5 dex, and potas- 
sium is under produced by about 2 dex. Scandium and 
titanium are under produced by about 0.3 dex. Cobalt 
is under produced by less than 0.5 dex. But for all other 
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Fig. 9. — Best-fit yields compared to observations of abundance patterns in the three most metal-poor stars. The abundances in 
HE0557-4840 are well-reproduced by the yields from model zl5G. Abundances of elements heavier than Na are modelled well by zl5D 
for HE0107-5240 and bracketed by models zl5D and z25D for HE1327-2326, but these models do not produce the high C,N, O, and Na 
abundances observed in the most metal-poor stars. Enrichment by a companion AGB star might explain these excesses. 

elements our synthetic yields fall within the error bars, or 
our yields from explosions of 1.2 or 2.4 Bethc bracket the 
data, demonstrating that SNe with energies between 1.2 
and 2.4 Bethe can account for the observed abundances. 
However, 15 M Q stars are needed in our IMF average 
to match the abundance patterns in both the Lai and 
Cayrel data. Spherical explosions with energies above 
2.4 Bethe are not required to match abundance patterns 
in the EMP data sets. While this does not rule out the 
existence of hypernovae, it does suggest that they are 
not required to explain the abundances in EMP stars. 
Also, as we noted earlier, our estimates of mixing are 
a lower bo und. With mor e detai led simulation physics 
like that in iKifonidis et al.l (|2006T >. it is likely that more 
heavy elements would be produced at a given explosion 
energy and that lower energy explosions alone would be 
sufficient to match the EMP data. 

Our IMF-averaged yields for rotating, zero-metallicity 
stars come closer to the observed abundances in EMP 
stars than do the one-d ime nsional parame t rizatio ns of 
iHeeer fe Wooslevi (|2008t ) or llwamoto et ail (|2005f ). the 
two-dimension a l non -rotating spherical exp losions of 
iJoggerst et al.l (|2009T > or the jet mo dels of ITominagal 
(12009)). Theangle-averaged yields of ITominagal (|2009f ) 
also underestimate potassium by nearly 2 dex. Nitro- 
gen, aluminum, scandium, and cobalt are also under 
produced, to a far greater degree than in our models. 
Although the jet-induced explosion may better account 
for the abundances in HE0107-5240 and HE1327-2326, 



an IMF average of spherical explosion yields of rotat- 
ing ze ro-m etallicity stars is a better fit to the lLai et al.l 
(|2008|) and lCavrel et al.l (j2004ft surveys. These two HMP 
stars may well be singular, rare objects. The more 
comprehensive surveys of EMP stars likely reflect more 
typical conditions and hence typical early supcrnovac. 
The good match of IMF- weighted averages of supernovae 
yields to these EMP star data suggests that supernovae 
of w 15 M with moderate explosion energies of less than 
2.4 Bethe were responsible for the bulk of early enrich- 
ment. Hypernovae may not be necessary to explain these 
abundance patterns. Our results imply that moderately 
sized, moderately energetic explosions produced most of 
the metals that went on to form these EMP stars. 

5.2. Incorporation of SNe Ejecta into Metal-Poor Stars 

Matching primordial nucleosynthetic yields directly to 
HMP and EMP abundances is problematic because not 
all of the elements may have been incorporated into 
subsequent metal-poor stars. Intervening hydrodynam- 
ical processes complicate the uptake of heavy elements 
into second-generation stars on several disparate spatial 
scales. First, as seen in our models, fallback onto the cen- 
tral remnant deep within the initial explosion can deprive 
later generations of certain elements, particularly Fe. On 
intermediate scales, HMP or EMP stars may have formed 
promptly from gas clumps created in hydrodynamical in- 
stabilit ies in primordial SN remnants (SNR) on 100 pc 
scales ([Whalen et al.l 1200 8) . If so, only the heavy ele- 
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Fig. 10. — IMF averages of yields compared to observations of EMP stars from lLai et al.l IMPOST ) and lCavrel et al.l J2004|). IMF averaged 
yields are performed over a single explosion energy (B, D, or G) and over all masses in the sample. Higher-explosion energy rotating Z = 
stars reproduce EMP abundances well. 15 Mq stars are needed to produce this good agreement with observations. 



mcnts that reach these instabilities by the time of their 
formation are those that become imprinted on the second 
stars. Their presence throughout the SNR can be traced 
back to their distribution in the early free expansion on 
scales of 10~ 3 - 10" 2 pc, which is well beyond those of 
our models. The velocity pro file of the free expansion i s 
nearly linear with radius (e.g. iTruelove fe McKedll999f ). 
so gas near the center of the explosion has much lower 
velocities than at greater radii. If this free expansion is 
not well mixed, heavier elements at small radii may never 
reach the formation sites of dense clumps at later times 
and appear in the resultant stars. Such differential mix- 
ing in the SNR can skew the abundance patterns of the 
second generation and is yet to be properly modeled. 

However, a rough idea of how these elements will be 
disseminated at later times can be obtained from their 
distribution with respect to velocity in our models, which 
we show in Figure 111! The velocities in this figure are 
those after all mixing has ceased, which occurs at a dif- 
ferent time for each star. Again, it is apparent that 



lower explosion energy SNe have greater fallback than 
do higher explosion energy SNe; their iron cores are es- 
pecially likely to fall back. For some stars, such as models 
u25B and z25B, u40B, u40D, z40B, and z40D, iron does 
not appear on the plot because it has all accreted onto 
the remnant. In general, because they exhibit more mix- 
ing, the z-series models show heavier elements mixed out 
to a higher percentage of the velocity of the lighter ele- 
ments than the u-series models. Nitrogen, because it is 
distributed throughout the helium envelope of the Z = 
models, is mixed out to a higher velocity than the other 
metals. 

Perhaps the most striking feature of these plots is that 
C, N, and O all consistently reach similar velocities that 
are noticeably greater than those of Fe. The separation 
between these elements and Fe will increase as the rem- 
nant expands homologously. This trend suggests that at 
later times regions of the ejecta at greater radii, which 
may be the sites of prompt second star formation, will be 
preferentially enriched with carbon, nitrogen, and oxy- 
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Fig. 11. — Distribution of elements with velocity. Progenitor mass increases from left to right; explosion energy increases towards the 
bottom of the figure; the first 3 rows correspond to Z = stars, the bottom 3 rows correspond to Z = 10~ 4 stars. The more compact 
Z = 10~ 4 Zq stars demonstrate more fallback (the long tails to the left of the peaks) than their Z = counterparts. Heavier elements are 
mixed our farther in velocity space for Z = than Z = 10~ 4 . Lower explosion energy leads to more fallback in all models, especially as 
mass increases. In all models in which Fe is ejected, it noticeably lags C, N, and O in velocity, implying that outer regions of the SNR will 
be preferentially enriched with those elements at later times. 

gen. Models zl5B, zl5D, and zl5G especially show C, 
N, and O separated in velocity from Mg and Fe. C, N, 
and O are also separated from Mg and Fe in mass coordi- 
nate, as can be seen in|6] This distribution may increase 
[CNO/Mg] while leaving [Fc/Mg] (as well as the ratios 
of element heavier than Mg to Mg) relatively unchanged, 
bringing our yields more into line with the abundances 
of HMP stars HE0107-5240 and HE1327-2326. Future 
simulations of the explosion over larger spatial scales and 
longer timescales are now being planned to evaluate these 
processes. 

6. CONCLUSIONS 

We have constructed two-dimensional simulations of 
the post-explosion hydrodynamics of rotating Population 
III and II core-collapse SNe. We find that the main ef- 
fect of rotation on the explosions is to alter the struc- 
ture of the stellar envelope through which they propa- 
gate. Even a small amount of rotation was enough to 
ex pand the compact st ationary Z = stars like those 
of iJoggerst et al.l (|2009f ) from blue to red. On the other 
hand, rotation had less effect on the presupernova struc- 
tures of Z = 10~ 4 Zq stars, which remained blue and dis- 
played a degree of mixing comparable to the non-rotating 
z-series models of IJoggerst et al.l (|2009f ). Once rotation 
was introduced in the progenitors, its actual magnitude 
had little bearing on the final structure of the star, or 
hence on the ejecta and yield of the SNe. There were no 



significant differences between R=5 and R=10 models of 
equal mass and explosion energy. 

The Z-series red supergiant models experience more 
mixing than the blue U-series models, and the redness or 
blueness of a supernova progenitor, not explosion energy, 
rotation rate, or mass, has the largest effect on stellar 
yields. We find that RT instabilities mix a greater per- 
centage of the star in the Z-series Population III models 
because the reverse shock takes longer to travel through 
the larger star, thus the prssure gradient is reversed for 
a longer period of time and the RT instability has more 
time to develop. The shorter travel time for the reverse 
shock through the more compact u-series models trun- 
cates the development of the RT instability in these stars 
after a shorter time than their z-series counterparts, re- 
sulting in less mixing. Likewise, there is less fallback onto 
the central remnant in the z-series explosions than in 
the blue u-series Population II models, whose envelopes 
are more tightly bound. At a given mass and metallic- 
ity, higher explosion energies lead to more mixing than 
lower energies, a trend that was most noticeable in the 
40 M Q runs. More massive stars exhibited less mixing 
than lower-mass stars. This was most apparent in the 
blue u-series stars, which experienced less mixing overall 
than their red z-series counterparts. Fallback increased 
with the mass of the progenitor and fell with increasing 
explosion energy. The 40 Mq models with explosion en- 
ergies of 1.2 Bethe or less did not eject any iron, nor did 
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25 M models with 0.6 Bcthe. 

Our rotating progenitor explosion models can account 
for the abundance patterns of the EMP stars and of 
one HMP star. The abundances in HE0557-4840, the 
least metal-poor of the HMP stars, with [Fe/H]=-4.75, 
are modeled well by the yields of a 15 M Q Z = SN 
with 2.4 Bethe. We are not able to replicate the high 
CNO to Fe ratios in the more metal-poor HMP stars 
HE0107-5240 and HE1327-2326, although we do repro- 
duce their elements above Na with Z = SNe mod- 
els between 15 and 25 M Q . Enrichment by an AGB 
companion may explain the extra C, N, O, and Na in 
HE 0107-5240 and HE1327-2326. Although they explain 
the abundance patterns in the two most iron-poor stars 
yet found, jet-driven explosions may not have been the 
most prevalent type of Population III supernova, just 
the one whose nucleosynthetic imprint survives in HMP 
stars in the Galactic halo today. The only stars with 
[Fe/H] < —5 that could have formed with masses low 
enough to be obse rvable today were those with high 
[C/Fe] and [O/Fc] dFrebel et all [20071 l200l . Another 
explanation for the high [CNO/Mg] ratios observed in 
these stars may lie in the final distribution of elements 
in our simulations. Models zl5D and z25D, provide the 
closest fits to elements heavier than Mg for HE0107-5240 
and HE1327-2326, respectively. These stars show C, N, 
and O separated from Mg and heavier elements in both 
mass and velocity coordinates. Preferential enrichment 
in these further-flung light elements may also reproduce 
the chemical abundance patterns in these stars without 
the need for a companion, whi ch has not been ob served 
despite continued observations (|Frebel et a l. 2008) . Low- 
mass stars with [Fe/H] < 5 might only have been able to 
form in a CNO-rich environment. 

The assertion that stars like the ones in this paper 
represent the dominant form of supernovae in the early 
Universe is bolstered by the fact that our yields are a 
better match to the EMP stars than HMP stars, which 
were likely enriched by a few generations of SNe whose 
cumulative C, N, and O imprint was sufficient to trigger 
a rollover from high-mass Population III stars to low- 
mass metal-poor Population II stars. Indeed, because the 
chemical makeup of EMP stars is the product of a few 
generations (and therefore a better cross section) of zero 
or very low metallicity explosions, the chemical abun- 
dances found in EMP stars may be a better metric of 
those of Z = and 10 -4 SNe in general. An IMF average 
over spherical-explosion SN yields from zero-metallicity 
rotating 15, 25, and 40 M Q progenitors with energies be- 
tween 1.2 and 2.4 B ethe repro d uces t he ab undance pat- 
terns found by the lLai et all (|2008[ ) and ICavrel et al.1 
(|2004ft survey of EMP stars with -4 <[Fe/H]< -3, with 
the exception of potassium, which is too low by ~ 2 dex 
in our models, and Al, which is too high by ~ 0.5 dex. 



The fact that 15 M progenitors must be included in 
the IMF average suggests that the metal-free stars that 
contributed the bulk of the metals to the early universe 
were of fairly low mass, extending down to the lower limit 
predicted for Population III stars. The abundance pat- 
terns of both HMP and EMP stars will soon come more 
fully into focus as the Sloan Extension for Galactic Un- 
derstanding and Exploration- II (SEGUE-II) component 
of the Sloan Digital Sky Survey-Ill (SDSS-III) uncovers 
hundreds of thousands of low-mass metal-poor stars in 
the Galactic halo. 

While the uncertainty about the explosion mechanism 
and geometry of core-collapse SNe means that these mix- 
ing estimates are likely lower bounds, our essential con- 
clusions are unchanged. Higher-mass stars may experi- 
ence more mixing, and less energetic blasts may eject 
amounts of heavy metals similar to those in our models 
with more realistic conditions for the first few seconds of 
the explosion. Hypernovae, or explosions with 10 times 
the energy observed in the nearby Universe, are not nec- 
essary to account for the abundances of the EMP stars: 
moderately enhanced explosion energies of less than 2.4 
Bethe suffice. Our models are consistent with observa- 
tions of abundances in EMP stars, suggesting that jet- 
driven mechanisms may not be necessary to explain these 
patterns. 

Our conclusions assume that SN ejecta that escapes 
fallback onto the remnant is uniformly mixed with the 
IGM and other SNR prior to incorporation into EMP 
stars, which may not be the case. Furthermore, even 
though the velocity distributions of elements in our zero- 
metallicity models are otherwise fairly uniform, Fe typ- 
ically lags C, N, and O. Metals in our u-series models 
exhibit more stratification in velocity space, suggesting 
that they will be differentially imprinted onto subsequent 
generations of stars. Simulations that follow the rem- 
nant out to greater spatial scales are necessary to test 
this hypothesis and are now under development. Mixing 
in three dimensions may also be different than in two, 
and will also be investigated in future simulations. 



Work at UCSC and LBL was supported in part by the 
SciDAC Program under contract DE-FC02-06ER41438. 
Work at LANL was carried out under the auspices of the 
National Nuclear Security Administration of the U.S. De- 
partment of Energy at Los Alamos National Laboratory 
under Contract No. DE-AC52-06NA25396. The simula- 
tions were performed on the open cluster Coyote at Los 
Alamos National Laboratory. Additional computing re- 
sources were provided on the Pleiades computer at UCSC 
under NSF Major Research Instrumentation award num- 
ber AST-0521566. 



REFERENCES 



Abel, T., Bryan, G. L., & Norman, M. L. 2000, ApJ, 540, 39 
— . 2002, Science, 295, 93 

Almgrcn, A. S., Beckner, V., Bell, J., Day, M., Howell, L., Joggcrst, 
C, Lijewski, M., Nonaka, A., Singer, M., & Zingale, M. 2010, 
ApJ 



Aoki, W., Frebel, A., Christlieb, N., Norris, J. E., Beers, T. C, 
Minezaki, T., Barklem, P. S., Honda, S., Takada-Hidai, M., 
Asplund, M., Ryan, S. G., Tsangarides, S., Eriksson, K., 
Steinhauer, A., Dcliyannis, C. P., Nomoto, K., Fujimoto, M. Y., 
Ando, H., Yoshii, Y., & Kajino, T. 2006, ApJ, 639, 897 

Beers, T. C. & Christlieb, N. 2005, ARA&A, 43, 531 

Bromm, V., Coppi, P. S., & Larson, R. B. 1999, ApJ, 527, L5 

— . 2002, ApJ, 564, 23 



17 



Burrows, A., Livne, E., Dessart, L., Ott, C. D., & Murphy, J. 2006, 
ApJ, 640, 878 

Calder, A. C, Fryxell, B., Plewa, T., Rosner, R., Dursi, L. J., Weirs, 
V. G., Dupont, T., Robey, H. F., Kane, J. O., Remington, B. A., 
Drake, R. P., Dimonte, G., Zingale, M., Timmes, F. X., Olson, 
K., Ricker, P., MacNeice, P., & Tufo, H. M. 2002, ApJS, 143, 
201 

Cayrel, R., Dcpagne, E., Spite, M., Hill, V., Spite, F., Francois, 
P., Plez, B., Beers, T., Primas, F., Andersen, J., Barbuy, B., 
Bonifacio, P., Molaro, P., & Nordstrom, B. 2004, A&A, 416, 
1117 

Christlieb, N., Gustafsson, B., Korn, A. J., Barklem, P. S., Beers, 
T. C, Bessell, M. S., Karlsson, T., & Mizuno-Wiedner, M. 2004, 
ApJ, 603, 708 

Ekstrom, S., Meynet, G., Chiappini, C, Hirschi, R., & Maeder, A. 
2008, A&A, 489, 685 

Frebel, A., Aoki, W., Christlieb, N., Ando, H., Asplund, M., 
Barklem, P. S., Beers, T. C, Eriksson, K., Fechner, C, Fujimoto, 
M. Y., Honda, S., Kajino, T., Minezaki, T., Nomoto, K., Norris, 
J. E., Ryan, S. G., Takada-Hidai, M., Tsangarides, S., & Yoshii, 
Y. 2005, Nature, 434, 871 

Frebel, A., Collet, R., Eriksson, K., Christlieb, N., & Aoki, W. 
2008, ApJ, 684, 588 

Frebel, A., Johnson, J. L., & Bromm, V. 2007, MNRAS, 380, L40 

— . 2009, MNRAS, 392, L50 

Fryxell, B., Arnett, D., & Mueller, E. 1991, ApJ, 367, 619 

Greif, T. H, Johnson, J. L., Bromm, V., & Klessen, R. S. 2007, 

ApJ, 670, 1 
Hamuy, M. 2003, ApJ, 582, 905 
Heger, A. & Woosley, S. E. 2002, ApJ, 567, 532 
— . 2008, ArXiv e-prints 

Herant, M. & Woosley, S. E. 1994, ApJ, 425, 814 
Hirschi, R. 2007, A&A, 461, 571 

Iwamoto, N., Umeda, H., Tominaga, N., Nomoto, K., & Maeda, K. 

2005, Science, 309, 451 
Joggerst, C. C, Woosley, S. E., & Heger, A. 2009, ApJ, 693, 1780 
Karlsson, T., Johnson, J. L., & Bromm, V. 2008, ApJ, 679, 6 
Kifonidis, K., Plewa, T., Janka, H.-T., & Muller, E. 2003, A&A, 

408, 621 



Kifonidis, K., Plewa, T., Scheck, L., Janka, H.-T., & Muller, E. 

2006, A&A, 453, 661 
Kitayama, T. & Yoshida, N. 2005, ApJ, 630, 675 
Lai, D. K., Bolte, M., Johnson, J. A., Lucatello, S., Heger, A., & 

Woosley, S. E. 2008, ApJ, 681, 1524 
Lohner, R. 1987, Comput. Methods Appl. Mech. Eng., 61, 323 
Mackey, J., Bromm, V., & Hernquist, L. 2003, ApJ, 586, 1 
Marek, A. & Janka, H.-T. 2009, ApJ, 694, 664 
Murphy, J. W. & Burrows, A. 2008, ApJ, 688, 1159 
Murphy, J. W., Ott, C. D., & Burrows, A. 2009, ArXiv e-prints 
Nakamura, F. & Umemura, M. 2001, ApJ, 548, 19 
Nomoto, K., Tominaga, N., Umeda, H, Kobayashi, C, & Maeda, 

K. 2006, Nuclear Physics A, 777, 424 
Norris, J. E., Christlieb, N., Korn, A. J., Eriksson, K., Bessell, 

M. S., Beers, T. C, Wisotzki, L., & Reimcrs, D. 2007, ApJ, 670, 

774 

O'Shea, B. W. & Norman, M. L. 2007, ApJ, 654, 66 
Scheck, L., Kifonidis, K., Janka, H.-T., & Muller, E. 2006, A&A, 
457, 963 

Tominaga, N. 2009, ApJ, 690, 526 

Tominaga, N., Umeda, H., & Nomoto, K. 2007, ApJ, 660, 516 
Truelove, J. K. & McKee, C. F. 1999, ApJS, 120, 299 
Umeda, H. & Nomoto, K. 2002, ApJ, 565, 385 
— . 2003, Nature, 422, 871 
— . 2005, ApJ, 619, 427 

Weaver, T. A., Zimmerman, G. B., & Woosley, S. E. 1978, ApJ, 
225, 1021 

Whalen, D., van Veelen, B., O'Shea, B. W., & Norman, M. L. 2008, 

ApJ, 682, 49 
Wise, J. H. & Abel, T. 2007, ApJ, 671, 1559 
Woosley, S. E. 1988, ApJ, 330, 218 

Woosley, S. E., Heger, A., & Weaver, T. A. 2002, Reviews of 

Modern Physics, 74, 1015 
Woosley, S. E. & Weaver, T. A. 1995, ApJS, 101, 181 
Yoshida, N., Oh, S. P., Kitayama, T., & Hernquist, L. 2007, ApJ, 

663, 687 



